Holstein magneto-polarons: from Landau levels to Hofstadter butterflies 
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We study the Holstein polaron in transverse magnetic field using non-perturbational methods. At 
strong fields and large coupling, we show that the polaron has a Hofstadter spectrum, however very 
distorted and of lower symmetry than that of a (heavier) bare particle. For weak magnetic fields, 
we identify non-perturbational behaviour of the Landau levels not previously known. 

PACS numbers: 71.38.-k,71.70.Di,73.43.Cd 
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The single polaron is the quintessential example of 
a dressed quasiparticle: as the electron interacts with 
bosonic modes from its environment, such as phonons, 
magnons or orbitons, it becomes "dressed" by a cloud 
of bosonic excitations. The properties of the resulting 
composite object - the polaron - can be significantly 
renormalized as compared to those of the bare parti- 
cle. Accurate numerical [l[ and analytical 0, Q ways 
to study such problems for any strength of the electron- 
boson coupling and in various dimensions have been de- 
veloped in recent years. This is to be contrasted with 
the case of dressed quasiparticles in strongly interacting 
systems, whose clouds consist of particle-hole excitations. 
Except for the few models with known exact solutions, 
their study away from perturbational regimes is still ham- 
pered by lack of accurate and efficient methods. 

Even though it is known that polarons have complex 
spectra, with substantial weight up to quite high energies 
above the low-energy polaron band, it is quite customary 
to expect that their behavior can be understood by think- 
ing of them as bare particles with a renormalized mass 
m*. In this Letter, we test this assumption by studying 
the response of polarons on two-dimensional (2D) lat- 
tices to an applied transverse magnetic field B. Note 
that for weak electron-phonon coupling, this problem has 
been studied extensively in continuous (as opposed to 
lattice) models using perturbation theory, because of its 
relevance to magnetotransport in 2D hetero-structurcs 
For weak B fields, it confirms the above-mentioned 
expectation by finding that the cyclotron frequency is de- 
fined by the polaron effective mass m* . As strong fields, 
it predicts an "un-dressing" of the quasiparticle and a 
cyclotron frequency controlled by the bare mass m. 

We use accurate non-perturbational methods to study 
the lattice problem for both weak and strong electron- 
boson coupling. To the best of our knowledge, this is 
the first time that a polaron lattice model in a trans- 
verse magnetic field has ever been investigated non- 
perturbationally. For weak coupling and weak fields we 
confirm the results of perturbational studies at low ener- 
gies. However, at higher energies we show the emergence 
of a complex pattern not predicted before, which is due 
to higher energy features of the polaron spectrum. 

For large magnetic fields and strong couplings, we in- 



vestigate for the first time the Hofstadter spectrum of 
small polarons. As is well known, if the flux (j) = Ba^ 
through the unit cell of a square lattice with lattice con- 
stant a is <p/4>o — p/q, where </>o = h/e is the quantum 
of magnetic flux and p and q are mutually prime inte- 
gers, the free-particle band splits into q sub-bands - the 
Hofstadter butterfly (5[ . We show that this splitting into 
q subbands holds for the small polaron as well. How- 
ever, the pattern is signiflcantly distorted and of lower 
symmetry than that of the bare particle, even for a cou- 
pling so large that m*/m « 100. This disagrees with 
the strong-coupling perturbational prediction of a simple 
mass renormalization. Taken together, these results show 
that at higher energies and / or for intermediary electron- 
boson couplings, the behavior of polarons is quite differ- 
ent from that of bare particles with larger mass m* . 

Model: We investigate the Holstein model @ - the sim- 
plest and most studied lattice model of electron-phonon 
interactions. The method we use is the Momentum Av- 
erage (MA) approximation, which has been shown to be 
highly accurate [3] not only for this but also for many 
other models, eg. with complex lattices, g{q) and g{k,q) 
coupling and disorder or inhomogeneities 0- Here 
we show that MA can also treat magnetic fields without 
any further approximations. The Hamiltonian is: 
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where i indexes sites on a square lattice and Ci, bi 
are electron/boson annihilation operators. The nearest- 



neighbour (nn) hopping ti 



exp 



has 



a Peierls phase defined by A{r) = ^{—y,x), VI is the en- 
ergy of the Einstein bosons and g is the strength of the 
electron-boson coupling. For B — Q, the spin of the elec- 
tron is irrelevant and is customarily ignored. For finite 
B, the spin degree of freedom is responsible for a trivial 
Zeeman splitting between spin-up and spin-down polaron 
states, which we also ignore in the following. 
The quantity of interest is the Green's function: 
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where |0) is the vacuum, G{lS) = [Tiw — "H + ^77] -'^ is the 
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resolvent, 77 > is infinitesimally small, and the second 
equality is the Lehmann representation in terms of the 
single-electron eigenstates 'H\a) — Ea\a). In particular, 
we will focus on the density of states (DOS): 

p{l^) = --Im G{i,i,oj) = V \{0\c^\a)\^S{nuJ-Ea). (2) 



Because this Hamiltonian is invariant to translations, 
there is no difference between local and total DOS. 

The MA approach has been discussed at length else- 
where @i S review here only the salient points. 
We start with the MA'"' formulation, which is equivalent 
to a variational expansion |a) = J2i j n't'i:3,-ncl{bj)'^\0) - 
i. e. SL cloud with any number of phonons can form at ar- 
bitrary distances from the electron, but all phonons are 
restricted to be at the same site [1]. Using this, we gener- 
ate equations of motion linking G{i,j,ui) to the general- 
ized Green's functions F„(i,j,a;) = (0|ciG(a;)c](&])"|0), 
as shown in Ref. The first (exact) equation reads: 

G(i, j, to) = Go(z, j, Lo)+gY,Fi ih h uj)Go{l, j, w). (3) 

I 

For any n > 1, we find within MA^^^ that Fn{i,j,uj) = 
gGo{j,j,uj-nn) [nFn-i{i,j,uj) + K+i(i,j, w)]. This re- 
currence equation is solved in terms of continuous frac- 
tions [3] to give Fn{i,j,Lj) — An{Lj)Fn-i{Lo), whcrc 
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G{i,j,uj) = Go{i,j,u~ Ej,f^(o)(w)) 



are independent of j because Tio is invariant to transla- 
tions. Using G{i,j,uj) — Fo{i,j,uj) in Eq. ^ gives: 

(5) 

where I]^^^(o) (lo) ~ gAi{uj). The only difference between 
this and the B = solution is that here Go{i,j,uj) is the 
free-electron propagator in the transverse magnetic jield. 
This can be calculated efficiently as shown in Ref. 9']. 

While MA(°) is accurate in describing ground-state 
(GS) properties for any effective coupling A = (7^/(4if2) 
so long as one avoids the extreme adiabatic limit ft/t ^ 

0, it does not properly account for the polaron-|-one- 
boson continuum that starts at Eqs + where Eqs is 
the polaron GS energy. This feature in the spectrum is 
due to excited states with a boson far away from the po- 
laron. To properly describe it, one needs to use or 
a higher level Q. At the MA*^^) level, the variational ba- 
sis is augmented with states like cj(6])"6; |0) with / ^ j, 

1. e. precisely the states contributing to the continuum. 
The equations of motions now also involve generalized 
Green's functions related to these states, which can be 
solved similarly like for the B = case Q. The final 
result is similar to Eq. ([5]), but the self-energy has the 
more accurate expression: 



(CO) 



5^Go(j, j,w) 



P(co) 

0.6 - 

0.4 - 
0.2 - 



(a) 



x5l 



u 



U 



-4.25 



-4.00 -3.75 

CO 



-3.50 



CO 

n 

-3.6 
-3.8 
-4.0 
-4.2 



1 1 1 1 


1 1 1 1 1 1 1 1 1 1 1 1 1 1 M>-l 






. < < 1 


(c) - 

. < < 1 < < < < 1 < < . < 1 < < . 



10 15 20 25 
n 



p(co) 

0.6- 

0.4 - 
0.2 - 
O.i 



(b) 



i 



x5 



)mm 




25 -4.00 -3.75 
CO 



-3.50 



0.00 0.02 0.04 0.06 0.08 0.10 



l~gGoUJ,^) [A^ito) ^ A,{ij - n)] 



FIG. 1; (color online) Density of states p{io) vs. energy u, 
for (a) A = and (b) A — 0.2. Dashed lines are for B — 
(times 5, for visibility), and full lines for 4'/4'o ~ 0.005. 
The arrow in (b) marks the edge of the polaron-|-one phonon 
continuum, (c) Energies uj„ of the Landau levels vs. n for 
(/>/(/>o ~ 0.0025, A = 0.2. The solid line is the prediction of 
Eq. ((B)); (d) The energy of the four lowest Landau levels as a 
function of (f>/(po for A = 0.2. The lines are the predictions of 
Eq. ([6]). Other parameters are t — 1,0, — 0.5, r] — 0.002. 



where uj — uj — CI — S^./^co) (w — 17). Again, the only 
difference from the B = result is the Go{j,j,uj) value. 
The self-energy's dependence only on uj is due to the 
simplicity of the Holstein model [2| . It becomes (weakly) 
non-local from the MA*^^-' level. Models with g{q) and 
gffc, q) coupling have strong momentum dependence in S 
[3| , but a finite B also only requires replacing free electron 
propagators with those in transverse field. 

All the results shown below are for the MA^^^ level. 
Like at i? = 0, this method is also equivalent to a summa- 
tion of all diagrams in the self-energy, up to exponentially 
small terms discarded from each. The resulting Green's 
function satisfies exactly the first 8 spectral weight sum 
rules, and with good accuracy the higher order ones [2j. 
While we do not know of any finite B numerical results 
for a direct comparison, the fact that the field is exactly 
included in the free propagator together with the argu- 
ments listed above, give us confidence that MA remains 
at least as accurate at finite i? as it is at B = 0. 

Results: We begin with a weak field and weak electron- 
boson coupling, where we can compare with known per- 
turbational results In Fig. [TJa) we plot the DOS 
with/without (full/dashed line) a very small field (/)/0o = 
0.005, in the absence of electron-boson coupling A = 0. 
The B = DOS is increased 5-fold for ease of view. As 
expected, it has a sharp rise at —4t and then increases 
slowly. For i? 7^ 0, we see the Landau levels as a succes- 
sion of Lorentzian peaks with width defined by rj. 
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Fig. shows the DOS at a weak couphng A = 0.2. 
The B = band-edge has moved below —At, due to the 
formation of the polaron band. The top of the polaron 
band and the jump marking the edge of the polaron-|-one- 
boson continuum at Eqs + ^ are clearly visible (arrow) . 
For B ^ 0, the polaron band splits into LLs with smaller 
spacing. Fig. [Ijc) shows their energies a;„ when (/)/4>o = 
0.0025. The line shows the perturbational prediction: 
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hujn = Egs + ^i"^* 
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where uj* = eB/m* is the cyclotron frequency, and we 
used the B — {) value of polaron effective mass, m* and 
of Egs- (Here, Eqs = -4.1288t and TO*/m = 1.0765). 
At low energies the agreement is very good, but it wors- 
ens as the LLs approach the continuum. Note that as 
expected, as the spacings decrease near the top of the po- 
laron band, so does the spectral weight in each LL. Even 
for low LLs, the agreement is worse at larger B, as shown 
in Fig. [l{d). The solid lines are Eq. ^ using again to*. 
The dashed lines show fits using a m*{B) = ■m*{l + jB), 
a correction to the polaron effective mass predicted by 
perturbation theory 4]. The agreement is much better 
although 7 increases with n, it is not a constant. Never- 
theless, we conclude that at low-energies the agreement 
with perturbation theory predictions is good [10|. 

At higher energies, however, it is not. Fig. [TJb) shows 
very different DOS in the continuum than below it. The 
failure of perturbation theory here is not surprising. Refs. 
Q use the free electron part as the large component, 
while electron-|-one-boson states are the small perturba- 
tion in the wavefunction. This is an accurate description 
at low energies but it fails at the top of the polaron band 
and inside the continuum, where the electron+one-boson 
states are dominant (for small A) while the free electron 
part is small. This failure of non-degenerate perturbation 
theory at these higher energies is well known for B — 



m. 



models, see for example Fig. 4 of Ref. 

What happens at these higher energies and also higher 
fields is shown in Fig. [2j the polaron LLs move to higher 
energies as B increases, until reaching an avoided crossing 
at an energy defined by the continuum band-edge as B ^ 
0, and which also moves higher with B. Above it, we see 
a whole sequence of such avoided crossings at energies 
that increase faster and faster with increasing B. 

The reason for this beautiful spectrum is easy to find. 
As mentioned, the B = polaron-|-one-boson continuum 
is due to excited states with a boson far from the polaron. 
At finite B, the continuum splits in a set of excited dis- 
crete states of energy fiujn + ^, each with a boson far from 
the polaron in a LL state. The DOS weights the spec- 
trum with the overlap with a free particle (zero bosons) 
state, see Eq. ([2]), so it vanishes at these energies. This 
explains the sequence of avoided crossings that occur at 
energies above that of the low-energy LLs. 
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FIG. 2: Contour plot of the density of states p{ij) vs lj and 
0/0O e [0.005, 0.06], for A = 0.2, t^l,Q^ 0.5, j? = 0.005. 



To our knowledge, this is the first time that this com- 
plex pattern is revealed. Perturbation theory Q predicts 
an avoided crossing at hujc ~ i.e. at (j^/fj^o ~ 0.04 for 
the values of Fig. [2] This is wrong but not surprising 
since as mentioned, non-degenerate perturbation theory 
is no longer valid at these energies. Perturbation the- 
ory also makes predictions about high fields huc 3> il. 
Here, Hofstadter butterfly effects become important for 
our model (they are absent in Refs. Q which use con- 
tinuous models). For small A, the DOS is quite complex 
because of overlap with the continuum and higher energy 
features. The results will be discussed elsewhere. 

Instead, here we focus on another interesting question, 
namely how like a particle is a strongly dressed quasipar- 
ticle? To answer this, we look at the Hofstadter spectrum 
of a small polaron, for A > 1. As is well-known, at i? = 
the small polaron band flattens considerably and a gap 
opens between it and the higher-energy features 
This gap allows us to look at the polaron response alone, 
avoiding overlap with these higher energy features. 

Fig. [3] shows results for A = 1.2, a valueiust above 
the crossover into the small polaron regime [J 0] . Panel 
(a) is the polaron band DOS at B = 0. The GS energy 
is signiflcantly lower because of the much larger binding 
energy, and the bandwidth is very narrow because of the 
large effective mass to*. As mentioned, this band is now 
separated by a gap from higher energy features. 

Panels (b)-(f) show the low-energy DOS for a magnetic 
flux (f>/(l)o = 1/6^1/2. The B = band indeed splits 
into q subbands for if)/4>o ~ p/q, as expected for a bare 
particle. For larger q the subbands become narrower and 
a smaller rj is needed. In panel (b) the thin line shows 
the DOS for rj — > 77/5. The peaks increase much less than 
5 times, proving that these are true continua (although 
very narrow and thus not yet fully converged at this 77), 
not discrete Lorentzians. For smaller q the subbands be- 
come much wider than rj and are already converged. 

We have checked (not shown) that, as required, spectra 
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FIG. 3: (color online) Density of states p{io) vs. energy uj, 
for A = 1.2, t = 1, Q — 0.5 and various magnetic fluxes 4>/4'o- 
The thick black lines correspond to 77 = 5 • 10~®, while the 
thin red line in panel (b) is for rj = 10""". The polaron band 
shows the Hofstadter signature, splitting into q sub-bands if 
t>Q — p/q. Higher energy features are not shown. 



are unchanged ii (/^/(po 1 ± '^/'/'o- We also see that 
Egs{B) increases significantly with (j)/4>o, reaching its 
maximum at (f>/(f>o — 1/2, consistent with the Hofstadter 
spectrum of a particle on a square lattice However, 
there are also big differences. The polaron spectra are 
asymmetric: the lowest subband has most of the weight 
and is quite distinct from the other subbands. This is 
very unlike the Hofstadter spectrum of the bare particle 
on a square lattice, which has particle-hole symmetry. 

This symmetry is lost even at i? = 0, where the DOS 
is not symmetric about the center of the band. The rea- 
son (see Fig. 2]) is that while the van-Hove singularity is 
still due to the flat i?(k) along the (0, tt) — (tt, 0) line, it 
is now located just below the upper band edge. More- 
over, as shown in panel (b), the quasiparticle weight is 
~ 2 orders of magnitude smaller here than near k = 0. 
Taken together, these explain the skewed shape of the 
B = DOS. They also show that nn hopping with a 
t* /t = m/m* = exp(— g^/ri^) as predicted by first order, 
strong-coupling perturbation theory [l2j, is not enough 
to fit i?(k), even though m*/m « 91. Second order per- 
turbation adds second and third nn hopping [12| : 
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which also give a poor fit, with non-monotonic behavior 
along all cuts shown in Fig. 21 except the (0, tt) — (tt, 0) 
line whose flatness is preserved. Indeed, to reasonably fit 
i?(k) one needs to add cos{nk,j;a + mkya) terms with up 
to |n| -|- |m| w 6. In other words, one needs to include 
terms at least up to 6th order in perturbation theory in 
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FIG. 4: Polaron energy E(k) (top) and quasiparticle weight 
Z{k) (bottom) vs k, for B = 0, t = 1, ^ = 0.5, A = 1.2. 



the hopping Hamiltonian to properly describe it. 

The long range hopping in E(k.) and the varying ^(k) 
explain the asymmetry of the polaron Hofstadter spec- 
tra. Taken together with the low-A results, they also show 
that a polaron is not behaving just like a bare heavier par- 
ticle with mass m*. Instead, its composite structure and 
the existence of higher energy states signal their existence 
in its finite- i? response. This has obvious implications for 
the interpretation of experimental data. 
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